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CN ■ ABSTRACT 

Observations of surface emission from isolated neutron stars (NSs) provide unique 

• challenges to theoretical modeling of radiative transfer in magnetized NS atmospheres. 

Recent work has demonstrated the critical role of vacuum polarization effects in de- 

fSJ ' termining NS spectra and polarization signals, in particular the conversion of photon 

, modes (due to the "vacuum resonance" between plasma and vacuum polarization) 

propagating in the density gradient of the NS atmosphere. Previous NS atmosphere 
models incorporated the mode conversion effect approximately, relying on transfer 

• equations for the photon modes. Such treatments are inadequate near the vacuum 
\^ I resonance, particularly for magnetic field strengths around B ^ Bi ~ 1 y. G, 

■ where the vacuum resonance occurs near the photosphere. In this paper, we provide 

an accurate treatment of the mode conversion effect in magnetized NS atmosphere 
' models, employing both the modal radiative transfer equations coupled with an ac- 

\ curate mode conversion probability at the vaccum resonance, and the full evolution 

equations for the photon Stokes parameters. In doing so, we are able to quantitatively 
calculate the effects of vacuum polarization on atmosphere structure, emission spectra 
and beam patterns, and polarizations for the entire range of magnetic field strengths, 
O I B — 10^^ — 10^^ G. In agreement with previous works, we find that for NSs with 

-4— > ' magnetic field strength B >2Bi, vacuum polarization reduces the widths of spectral 

I features, and softens the hard spectral tail typical of magnetized atmosphere models. 

' For B < Bi/2, vacuum polarization does not change the emission spectra, but can 

. \ significantly affect the polarization signals. Our new, accurate treatment of vacuum 

' polarization is particularly important for quantitative modeling of NS atmospheres 

H \ with "intermediate" magnetic fields, B ~ (0.5 — 2)B;. We provide fitting formulae 

- - ■ for the temperature profiles for a suite of NS atmosphere models with different field 

strengths, effective temperatures and chemical compositions (ionized H or He). These 
analytical profiles are useful for direct modeling of various observed properties of NS 
surface emission. As an example, we calculate the observed intensity and polarization 
lightcurves from a rotating NS hotspot, taking into account the evolution of pho- 
ton polarization in the magnetosphere. We show that vacuum polarization induces a 
unique energy-dependent linear polarization signature, and that circular polarization 
can be generated in the magnetosphere of rapidly rotating NSs. We discuss the im- 
plications of our results to recent observations of thermally emitting isolated NSs and 
magnetars, as well as the prospects of future spectral and polarization observations. 

Key words: magnetic fields - radiative transfer - stars: atmospheres - stars: magnetic 
fields - stars: neutron - X-rays: stars. 
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1 INTRODUCTION 

Thermal surface emission from isolated neutron stars (NSs) (e.g., iKaspi et allboosll can potentially provide invaluable infor- 
mation on the physical properties and evolution of NS equations of state at super-nuclear densities, cooling histories, magnetic 
fields, and surface compositions (see, e.g., Prakash et al. 2001; Yakovlev & Pethick 2004 for review). In recent years, significant 
progress has been made in detecting such radiation using X-ray telescopes such as Chandra and XMM-Newton. For example, 
the spectra of a number of radio pulsars (PSR B1055-52, B0656-I-14, Geminga and Vela) have been observed to possess thermal 
components that can be attributed to emission from NS surfaces and/or heated polar caps (see, e.g., Becker & Pavlov 2002). 
Phase-resolved spectroscopic observations have become possible, revealing the surface magnetic field geometry and emission 
radius of the pulsar (e.g., Caraveo et al. 2004; De Luca et al. 2005; Jackson & Halpern 2005). Chandra has also uncovered 
a number of compact sources in supernova remnants with spectra consistent with thermal emission from NSs (see Pavlov et 
al. 2003), and useful constraints on NS cooling physics have been obtained (Slane et al. 2002; Yakovlev & Pethick 2004). 

Surface X-ray emission has also been detected from a number of soft gamma-ray repeaters (SGRs) and anomalous X-ray 
pulsars ( AXPs) — these are thought to be magnetars, which emit radiation powered by the decay of superstrong {B ^ 10^* G) 
magnetic fields (see Thompson & Duncan 1995,1996; Woods & Thompson 2005). The quiescent magnetar emission consists of 
a blackbody-like component with temperature (3-7) x 10® K, and a power-law tail from 2—10 keV with photon index 2 — 3.5 
(e.g., Juett et al 2002; Kulkarni et al. 2003; Tiengo et al. 2005), as well as significant fiux at ~ 100 keV (Kuiper et al. 2004; 
Kuiper et al. 2006). Somewhat surprisingly, the observed thermal emission does not show any of the spectral features one 
might expect, such as the ion cyclotron line around 1 keV (for typical magnetar field strengths). 

The seven isolated, radio-quiet NSs (the so-called "dim isolated NSs"; see Haberl 2005) are also of great interest. These 
NSs share the property that their spectra appear to be entirely thermal, indicating that the emission arises directly from 
the NS atmospheres, uncontaminated by magnetospheric processes. The true nature of these sources is unclear at present: 
they could be young cooling NSs, or NSs kept hot by accretion from the ISM, or magnetars and their descendants (see 
van Kerkwijk & Kulkarni 2001, Mori & Ruderman 2003; Haberl 2005; Kaspi et al. 2005). While the brightest of these, RX 
J1856. 5-3754, has a featureless spectrum remarkably well described by a blackbody (Drake et al. 2002; Burwitz et al. 2003), 
absorption lines/features at E ~ 0.2-2 keV have recently been detected from at least four sources, including IE 1207.4- 
5209 (at 0.7 and 1.4 keV, possibly also 2.1, 2.8 keV; Sanwal et al. 2002; Bignami et al. 2003; DeLuca et al. 2004; Mori et 
al. 2005), RX J1308.6-f 2127 (0.2-0.3 keV; Haberl et al. 2003), RX J1605. 34-3249 (0.45 keV : van Kerkwijk et al . 2004), RX 
J0720 .4— 3125 (0.27 keV; Haberl et al. 2004a), and possibly several additional sources (see lHaberl et aljEo04bl : IZane et alJ 
I2OO5I) . The identification of these features remains uncertain, with suggestions ranging from electron/ion cyclotron lines to 
atomic transitions of H, He or mid-Z atoms in a strong magnetic field (Ho & Lai 2004; Pavlov & Bezchastnov 2005; Mori et 
al. 2005). These sources also have different X-ray lightcurves: for example, RX J1856. 5-3754 and RX J1605-I-3249 show no 
variability (pulse fraction ^ 1-3%), while RX J0720-3125 shows single-peaked pulsations (P = 3.39 s) of amplitude ~ 11%, 
with the spectral hardness and line width varying with the pulse phase, and RX J1308-I-2127 shows double-peaked pulsations 
(P = 10.3 s) with amplitude ~ 18%. Another puzzle concerns the optical emission: for at least four of these sources, optical 
counterparts have been identified, however, t he optical flux is larger (by a facto r of 4-10) than the extrapolation from the 
blackbody fit to the X-ray spectrum (see, e.g.. iTrumoer et ahllioM: iHaberiboO.'il^ . 

The spectrum of NS thermal radiation is formed in the atmosphere layer (with scale height ^ 0.1 — 10 cm and density 
10"^ - 10^ g cm ) that covers the stellar surface. Thus, to properly interpret observations of NS surface emission, detailed 
modeling of NS atmospheres in strong magnetic fields is required. 



1.1 Previous work on magnetic neutron star atmosphere models 

The first magnetic NS atmos phere models were constructed by Shibanov et al. (1992) (see also lPavlov et al.ll995l : lRaiagopal et alJ 
ll997l:IZavlin fc Pavlovll2002l) . who focused on moderate field strengths B ~ 10^^-10^^ G and assumed full ionization (see also 
Zane et al. 2000 for atmosphere models with accretion).^ Similar ionized models for the magnetar field regime {B Z 10" G) 
were studied by Zane et al. (2001), Ozel (2001) and Ho & Lai (2001,2003). An inaccurate treatment of the free-free opacities in 
the earlier models (Pavlov et al. 1995) was corrected by Potekhin & Chabrier (2003), and the correction has been incoporated 
into later models (Ho et al. 2003,2004). Recent works (Lai & Ho 2002,2003a; Ho & Lai 2003) have shown that in the magnetar 
field regime, the effect of strong-field quantum electrodynamics significantly influences the emergent atmosphere spectrum. 
In particular, vacuum polarization gives rise to a resonance phenomenon, in which photons can convert from the high-opacity 
mode to the low-opacity one and vice versa. This vacuum resonance tends to soften the hard spectral tail due to the non- 
greyness of the atmosphere and suppress the width of absoprtion lines (see Lai & Ho 2003a for a qualitative explanation). 

An earlier attempt bv lMilleil 

^^23) adopted a polarization-averaging procedure for the radiative transport which is rather inaccurate 
since most of the photon flux is carried by the low-opacity photon mode. 
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Even for modest field strengths {B ^ 10^'* G), vacuum polarization can still leave a unique imprint on the X-ray polarization 
signal (Lai & Ho 2003b). ^ 

Because a strong magnetic field greatly increases the binding energies of atoms, molecules and other bound species (see Lai 
2001), these bound states may have appreciable abundanc es in the NS atmospher e (Lai & Salpeter 1997; Potekhin et al. 1999). 
Early considerations of partially ionized atmospheres (e.g.. lR,aiagopal et alll99^ relied on oversimplified treatments of atomic 
physics (e.g., ionization equilibrium and equation of state) and nonideal plasma effects in strong magnetic fields. Recently, a 
thermodynamically consistent equation of state and opacities for magnetized {B = 10^^ — lO^"" G), partially ionized H plasma 
have been obtained (Potekhin et al. 1999; Potekhin & Chabrier 2003,2004), and the effect of bound atoms on the dielectric 
tensor of the plasma has also be studied (Potekhin et al. 2004). These improvements have been incorporated into partially 
ionized, magnetic NS atmosphere models (Ho et al. 2003; Potekhin et al. 2004). Finally, for sufficiently low te mperatures 
and /or stro ng magnetic fi elds, the NS atmosphere ma y undergo a phase transiti on into a cond ensed state (see iRudermarJ 
ll97ll:ljoneJ.1986 : Neuhauser et ai]|l987l: Lai fc Salpeter,. 1997.: , Lai. D. .2001. : -Medin fc Laill2006allbl) . Thermal emission from 
such a condensed surface has been studied by van Adelsberg et al. (2005) (cf. TuroUa et al. 2004; Perez- Azorin et al. 2005). 



1.2 This paper 

While previous works have identified the importance and trend of the vacuum polarization effect (Lai & Ho 2002,2003a), the 
implementation of the effect in NS atmosphere models (Ho & Lai 2003,2004) has been based on approximations (see §2). 
So far all studies of magnetic NS atmospheres have relied on solving the transfer equations for the specific intensities of the 
two photon modes. As discussed in Lai & Ho (2003a) and reviewed in §2 below, these equations cannot properly handle the 
vacuum-induced mode conversion phenomenon because mode conversion intrinsically involves interference between different 
modes. 

In this paper we provide a new, quantitatively accurate treatment of vacuum polarization effects in radiation transfer 
for fully ionized NS atmospheres. Our work confirms the semi-quantitative results obtained in previous works based on an 
approximate treatment of vacuum resonance (Ho & Lai 2003; Lai & Ho 2003b). Moreover, our new treatment allows us to 
quantitatively predict the spectral and polarization properties of NS atmospheres with field strengths varying from 10^^ G to 
10^^ G. 

The remainder of our paper is organized as follows: 32] describes the basic physics of vacuum polarization in radiative 
transfer and resonant mode conversion; ^details two methods of solving the radiative transfer problem and our basic physics 
inputs, as well as tests of our numerical method; ^presents the basic results from our models, including atmosphere structures 
(with fitting formulae for the temperature profiles), spectra, and emission beam patterns; ^considers the observed polarization 
signals of atmosphere emission for a simple geometry (i.e., a rotating NS hotspot); and discusses the implications of our 
results. 



2 EFFECT OF VACUUM POLARIZATION ON RADIATIVE TRANSFER 

Before describing our quantitative treatment of the vacuum polarization effect in NS atmospheres, it is useful to summarize 
the basic physics of the effect (see also Lai & Ho 2003a) and discuss the limitations of previous treatments. 

Photons (with energy E Ebb = heB/nieC, the electron cyclotron energy) in magnetized NS atmospheres usually 
propagate in two distinct polarization states, the ordinary mode (denoted by "O") and the extraordinary mode (denoted by 
"X"), which are polarized (almost) parallel and perpendicular to the plane made by the magnetic field and direction of photon 
propagation, respectively. In strong magnetic fields, the dielectric tensor describing the atmospheric plasma of a NS must be 
corrected for QED vacuum effects (Gnedin et al. 1978; Meszaros & Ventura 1979; Pavlov & Shibanov 1979; Meszaros 1992). 
For a photon propagating in a medium of constant density p, the plasma and vacuum contributions to the dielectric tensor 
"cancel" each other out at a particular energy given by 

Ev = 1.02 (y; pi)'/" ^rZ/s keV, (1) 

where Ye — Z/A {Z, A are the atomic number and mass number, respectively), pi — p/{l g cm~"^), B14 — 5/(10^* G), and 
/s ~ 1 is a slowly varying function of B [see eq. (2.41) of Ho & Lai (2003)]. At the resonance, both modes become circularly 
polarized. A number of previous papers (e.g., Meszaros 1992) emphasized the sharp X-mode opacity feature associated with 

2 iBuIik fc Miiieil Jl997f) studied the effect of vacuum polarization on Compton scattering in an isotliermal magnetized plasma, in tlie 
context of SGR bursts. They did not calculate self-consistent atmosphere models and did not find the effects (suppression of spectral 
lines and hard tails) discussed in the works of Ho fc Lai. . Zanc e t al. (2001) included vacuum polarization effects in their models, but 
did not report any of the important effects of QED. IOze]l '" i2n0lih also included vacuum polarization effects, but came to the opposite 
conclusion (i.e., that vacuum polarization effects cause thermal spectral tails to become harder), which is incorrect. 
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Figure 1. Photon free-free absorption opacities for X and O polarization modes as a function of energy at i? = 10^* G, T = 10® K, 
SkB = ""/^ and p = 5.4 g cm~^. Vacuum polarization induces the sharp resonance feature for the X-mode opacity at Ey. This "spiky 
opacity" can affect the emergent radiation spectrum from magnetized NSs, but does not include all the effects associated with the vacuum 
resonance. 



the resonance (see Fig. . It might seem that to understand the vacuum polarization effect in radiative transfer, all one needs 
to do is to include this spike in the opacity (e.g., Ozel 2003). However, this is not the whole story. 

A more useful way to understand the effects of the vacuum resonance is to consider a photon with given energy E, 
traversing the density gradient of a NS atmosphere. The photon will encounter the vacuum resonance at the density 

pv = 0.96 Y-'e! bI^ fg^ g cm-^ (2) 

where Ei = E/{1 keV). lLai fc 53 |20o3) showed that the photon undergoes resonant mode conversion when the adiabatic 
condition E > Ead is satisfied, with' 



3 



Eaci = 2.h2[fBta.nekB\l-{EB^/Ef\f'^ , (3) 

where 6kB is the angle between the magnetic field and direction of propagation, Ebi ~ 0.63{Z/A) keV is the ion cyclotron 
energy, and Hp = \ds/d\np\ is the density scale height along the ray. Thus, an adiabatic 0-mode photon encountering the 
vacuum resonance will convert into an X-mode photon, and vice- versa. In general, for intermediate energies E ~ Ead, photons 
undergo partial conversion, in which an O-mode converts to a X-mode (and vice-versa) with probability 1 — Pjump, where 
Pjump is the non-adiabatic jump probability 



-fjump — exp 



liE/Eadf 



(4) 



Due to free-free absorption, the X-mode opacity is suppressed relative to the O-mode by a factor of {Ebc/ E)^ , where the 
electron cyclotron energy is Ebe = 1158 B14 keV; thus, the mixing of photon modes at the resonance can have a drastic effect 
on the radiative transfer. For magnetic field strengths satisfying (Lai & Ho 2003a, Ho & Lai 2004) 

B ^ B; ~ 6.6 X 10'^ TfT^/^E'^^'^S-^^'^ G, (5) 

where Tg = r/(10'' K) and S = l-e"'^''''^, the vacuum resonance density lies between the X-mode and O-mode photospheres 
for typical photon energies, leading to suppression of spectral features and softening of the hard X-ray tail characteristic of 
ionized hydrogen atmospheres. For "normal" magnetic fields, B Bi, the vacuum resonance lies outside both photospheres, 
and the emission spectrum is unaltered by the vacuum resonance, although the observed polarization signals are still affected 
(Lai & Ho 2003b). 

In their implementation of the vacuum resonance effect in NS atmosphere models. Ho & Lai (2003) considered two 
Since E^d depends on E, one needs to solve for E > E^d to determine the adiabatic region. See Fig. 6 of lLai fc Hoi |20033). 



© 2006 RAS, MNRAS OOO.ITllsTI 



Atmosphere Models of Magnetized Neutron Stars 5 



limiting cases: (i) complete mode conversion (Pjump ~ 0), which is equivalent to assuming that E 2> Ead is satisfied for all 
photon energies; (ii) no conversion (Pjump = 1), which is equivalent to assuming E <C Ead for all photons. In the former case, 
all X-mode photons are converted to the 0-mode at the resonance (and vice- versa), whereas in the latter, such conversion 
is neglected. In both cases, radiative transfer equations based on photon modes can be used, as long as one properly defines 
the modes across the resonance (Ho & Lai 2003). We expect that the complete and no conversion limits bracket the correct 
solution. In case (ii), one only has the "narrow spiky opacity" effect associated with the resonance. Lai & Ho (2002) estimated 
the width of this opacity spike and emphasized the importance of resolving the spike. In both limits, vacuum resonance has 
qualitatively the same effects on the emergent spectrum, i.e., suppression of lines and softening of hard spectral tails (Lai & 
Ho 2002; Ho & Lai 2003), although for B ~ (a few x 10^^) — 10^* G, appreciable quantitative differences in the spectra using 
the two limits are produced (Ho & Lai 2004). 

As mentioned before, all studies of radiative transfer in magnetized NS atmospheres so far have relied on solving the 
transfer equations for the specific intensities of the two photon modes (e.g., Meszaros 1992; Zavlin & Pavlov 2002). These 
equations cannot properly handle the vacuum-induced mode conversion phenomenon. In particular, photons with energies 
0.3-2 keV are only partially converted across the vacuum resonance (this is the energy range in which the bulk of the radiation 
emerges and spectral lines are expected for B ~ 10^* G). In addition, the phenomenon of mode collapse (when the X and 
0-modes become degenerate) occurs when dissipative effects are included in the plasma dielectric tensor, and the concomitant 
breakdown of the Faraday depolarization condit ion near the resonance fu rther complicates the standard treatment of radiative 
transfer based on normal modes. As shown bv lGnedin fc Pavlovl (^T^), the modal description of radiative transfer is valid 
only in the limit |Re(nx — wo) I ^ |Im(wx + no)\, where nx and no are the indices of refraction corresponding to the X 
and O-modes, respectively. IHo fc Lai (|20fl 3^1 showed that, for a narrow range of energies around the vacuum resonance, this 
condition can be violated, and the violation becomes especially pronounced in the magnetar field regime. It is not obvious 
whether the mode collapse significantly alters the radiative transfer. Thus, to account for the vacuum resonance effect in a 
quantitative manner, one must solve the transfer equations in terms of the photon intensity matrix (Lai & Ho 2003a) and 
properly take into account the probability of mode conversion. This is one of the main goals of our paper. 



3 METHOD 

3.1 Partial mode conversion using photon mode equations 

3.1.1 Radiative transfer equation 

For our models, we consider plane-parallel, fully ionized H or He atmospheres, with the magnetic field oriented normal to the 
surface. The standard method used in all previous work involves solving the coupled radiative transfer equations for the two 
modes of photon propagation. These are given by 

± = ^ [Hir, ±m) - Si{r, ±m)] (6) 

where /^(r, ji) is the specific intensity for mode j, ^ = k • z ^ 0, — + is the total opacity (with contributions from 
free-free absorption and scattering, see below), kt ~ 0.4 cm^ g~^ is the Thomson scattering opacity, r is the Thomson optical 
depth (defined by dr — —pKrdz), and Si is the source function, defined below. Eq. Q is solved subject to the constraints 
of hydrostatic and radiative equilibria, as well as constant radiative fiux Prad, given by: 

(7) 



HiT 

dv 
i=i 



^,f=(^_j^)=0, (8) 



2 

Prad = 2^y" / dv i d/i/i [li(pi) - lii^ii)] = a^^T^s (9) 



where P is the pressure of electrons and ions, g = (^^) (l — ^^^) ^''^ = 2.4 x 10^* cm s ^ is the surface gravitational 
acceleration (we adopt M = 1.4 Mq and ii = 10 km throughout the paper), Jl = (1/2) f^ d/i [Hi/J-) + /,:i(— /i)] is the mean 
specific intensity, Bi, is the Planck function, and Tcff is the effective temperature of the atmosphere. To integrate eq. @ subject 
to the conditions of f^-l j^. we assume the ideal gas equation of state for both protons and electrons; electron degeneracy 
effects are neglected. IHo fc Lail l)200lll showed that the effect of this approximation on the atmosphere is negligible. Note 
that in general, thermal conduction due to electrons also contributes to the total flux. However, we have found that in the 
atmosphere region of interest, the conduction flux is always less than a few percent of the total flux, and is therefore neglected. 
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3.1.2 Photon modes and opacities 

The properties of magnetized atmospheric plasma can be described by a complex dielectric tensor dCinzbur j|l964) . In a 
coord inate system with the magnetic field aligned with the z-axis, the plasma contribution to the dielectric tensor takes the 
form l)Lai fc Holl2003ah : 



.4 




[e^^'^] =\ -ig 6 , (10) 



where 

e±g ^ 1 ^ — —J — ^ — ^ — (11) 

(l+i7re±We )(l+i7"TMi ) + «7ei ) 

V ~ 1 H « • (12) 

l+«(7e.+7re) I + + Iri) 

In eqs. Illlll - 112| l we have defined the dimensionless ratios Ue = (EBe/E)^, Ui = {EBi/E)"^, Ve = (Epe/E)^, Vi = (Epi/E)^, 
where Epe = h{4:nnee^ /rrie)^^^ = 0.02871(ye pi)^^^ keV is the electron plasma energy, and Epi = {Zme / Amp)Epe = 6.70 x 
10"^ Ye p\l'^ keV is the ion plasma energy. The dimensionless damping rates 7^'" — v^^^ juj (for electron-ion coUisional 
damping), 7^^ = Vre/oj (for electron radiative damping), and 7^1 = Vrijio (for ion radiative damping) are given by, 

7i'" = 9-2xl0-^||^(l-e-/''«-),l,„ (13) 

7.e = 9.5 X lO^^Si, (14) 

-yri = 5.2 X lO"''-:-^;!. (15) 

The quantities and are the velocity-averaged magnetic Ga unt factors perpe ndicular and parallel to the magnetic field, 
respectively; they are calculated using eqs. (4.4.9)-(4.4.12) from iM'eszaro j l)l992l) .^ This calculation includes contributions 
fro m electrons in the gro und Landau level only. Gaunt factors including contributions from excited states have been derived 
bv lPotekhin fc Ghabrieil (Eo04il . Nevertheless, for energies well below Ebe, the differences between the two calculations are 
neghgible JPotekhinlboog) . 

Vacuum contributions to the dielectric tensor can be taken into account by making the following substitutions into the 
tensor of eq. llOU : 



e^e' = e + a— 1, r^^-q'^rj + a + q— 1, (16) 

where a and q are vacuum parameters given by the expressions in, e.g., Heyl & Hernquist (1997) and Potekhin et al. (2004) 
(the latter also contains general fitting formulae). Solving Maxwell's equations for the anisotropic medium yields two modes 
of propagation. In a coordinate system where the wave vector k is along the 2-axis and the magnetic field lies in the xz plane 
(such that k X B = — sin6'fe_By), the mode eigenvectors can be written as 

«± = (rT^^Ti^:7^(^^±'i'^^^±)' (17) 

where the ellipticity K± = —ie^/ey of mode ± is given by 



K±=p± VP^, (18) 

with r = 1-1- (m/a) sin'^ 9kB {m is another vacuum polarization parameter; Heyl & Hernquist 1997; Potekhin et al. 2004), and 
the polarization parameter [3 is 

o e'"-ff'-eV(l-|-m/a)sin"6'fcB 

^ 2gn' cose.B- ^ ' 

The ^-components of the mode eigenvectors are given by 

^ _ (e - - g) singfcB cos 6kBK± + gsinOkB ^2o) 
e sin^ 9kB + iv + q) cos^ OtB + a- 1 
Note that when the modes are labelled according to eq. IjlH^ . the K± vary continuously across the vacuum resonance 



Note that eq. (13) of lLai fc Hd i2003a^ should be 7^ = 7^(1-1- ) ~ y^^. This substitution should be applied to all the appropriate 
equations in lLai fc Hd J2003a|^^__^^ 

^ Note that eq. (4.4.12) oflMfezlTOi EUl) should be a± = (p ± [p'^ + 2mhui]'^^^)^{2mkBT)-''- . 
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(/3 — 0), and do not cross each other in the absence of dissipation jLai fc Holl2003al) . Another way of labeling the modes, 
commonly adopted in the literature (e.g. lMeszaroj|^^92^ . is 



i + i-iy 1 + 



r 



1/2 



(21) 



According to this labeling scheme, j = 1 corresponds to the X-mode {\Ki\ < 1) and j ~ 2 corresponds to the 0-mode 
{\K2\ > 1). Obviously, Ki and K2 are not continuous functions across the vacuum resonance. It is also clear that a given + 
mode (or - mode) which manifests as the X-mode (0-mode) before the resonance switches character after the resonance. 

Using the mode eigenvectors and the componenets of the dielectric tensor, expressions for the free-free absorption and 
scattering opacities can be obtained. The cyclic components of the mode eigenvectors in a rotating frame with the magnetic 
field along the z-axis are: 



(ex 



1 

Kj sin 9kB 



1 ± \Kj cos6kB + K^i sin 6*^31 



2(1 

|2 



1 + \KA^ + \K. 



(22) 



(23) 



Note that in the above expression, j indicates the mode, and the ± subscri pt should not be c onfused with the K± labeling of 
photon modes. The free-free absorption opacity for mode j can be written 



K 



with 

K± 

A± 

7± 



' A -L 
-i;eA±7ei, 



cp 

UJ I 



{l±uj ) [I^U/ ) +7± 
7^ + (1 ± uJ'^Yfri + (1 T u]'^)'yre. 



(24) 

(25) 
(26) 

(27) 
(28) 



Note that these expressions include the contribution of electron-ion Coulomb collisi ons to the free-free absorption opac ity in 
a consistent way. They correct the free-free opacity adopted in earlier papers (e.g. IPavlov et al."l995': 'Ho & Lai"200l') , and 
they agree with the correct expressions given in »Potekhin fc Chabrier. 12003.). an d th ose used by Ho et al.. (|2flfl3t . ■2004 ') . 



. . „ „ „.. py tip g 

The scattering opacity from mode j into mode i is given bv lventurJ |l979ll fsee lHo fc Laill200lt l: 



1 

Kfi = YeKT ^ [(1 + 



1/2x2 



'a 



EI- 



l/2s2 , 2 



Al 



(29) 



where je = Jei+lre, 7i ~ 7S+7ri, and Al, — (3/4) J^_^ dii\e]^\'^ . The total scattering opacity from mode j is then Kj'^ = k^^. 



3.1.3 Source function 

The source function in eq. @ can be written as 



4 



+ 



where iVentur Jl979l) 



2-K 



2 



[ii{p) + ni-p)] 



dn' 



8tv 



1 



+ 



Stt V Arrir, 



~A 



1 

1 - ccuV^ 



(30) 



(31) 



Following IHo fc LaJ i200ll) . it is a good approximation to assume that the differential scattering cross-section is independent 
of the initial photon direction. The resulting approximate source function is: 

4 (m) , Kfiin) cut 



SIM 



«:*°'(At) 2 



(m) 47r 



(32) 



where uj = (27r/c) J^_^ dp! ll{p!) is the specific energy density of mode j. Note that the source function depends on the specific 
intensity in all directions, and thus depends on the solution to the radiative transfer equation. Therefore, we calculate Si 
iteratively, according to the scheme described in H3.1.4tH3.1.5l 
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3.1.4 Solution to transfer equation for photon modes including partial mode conversion 

We describe above how vacuum polarization effects can be incorporated into the free-free absorption and scattering opacities 
for the phot on modes. Howev er, these opacity effects do not capture the essence of the vacuum resonance phenomena. As 
discussed in iLai fc Hoi J2003ah . solving the transfer eq. @ using K± [eq. III8I1 I as the basis for the photon modes amounts 
to assuming complete mode conversion (-fjump — 0), while using lKi^2 [eq. 1121^ ] corresponds to assuming no mode conversion 
(fjump = 1). This was the strategy adopted by previous works. To correctly account for the vacuum resonance effect, it is 
necessary to use the jump probability Pjump [eq. (jlj] to convert the mode intensities across the resonance according to the 
formulae: 

Ix fjump + (1 — Pjump )/o, (33) 

lo ^ Pjump /0 + (l-Pjump)/x. (34) 

Note that since the resonance density depends on photon energy, the standard Feautrier procedure for integrating the radiative 
transfer equation cannot be used here, as there is no simple way to incoporate eqs. H33|I - I|34|P into the method of forward and 
backward substitution employed by Feautrier (see Mihalas 1978, §6-3). Instead, we use the standard Runga-Kutta method to 
integrate the transfer eq. @ in the upward and downward directions starting from the boundary conditions: 

II{t r^ax,+y) (35) 

Ii{T ^ Train,— ^ (36) 

The Runga-Kutta integration is stopped at the resonance, where eqs. H33|l and 1134^ are used to convert the mode intensities. 
Then the integration is continued to completion. The limits r,nax and Tmin are chosen to span 5-8 orders of magnitude. 
This insures that (1) photons begin their evolution at densities greater than the X-mode decoupling depth, and (2) both 
X-mode and 0-mode photons are fully decoupled from the matter at the outermost edge of the atmosphere (r = Tmin). We 
finite-difference eq. Q as: 

±M^^^[«(/-5)+«'(/'-5')] (37) 
where / = I^t, ±fi), J' = I^t', ±^), k = Hif \ n' = 4°*(r'), S = Si{T, ±/i), S' = SUt' , ±^l), with r' = r + Ar. This gives 



1 T 



(38) 

\ /.llK^r I /.ll.Ksr 

This formula yields stable integrations whose results are not strongly dependent on grid spacing (see >13.41 . 

To summarize, our method for integrating the radiative transfer with partial mode conversion is as follows: (1) For given 
E and 6kB, we integrate eq. Q using 1381 from Tmax to the vacuum resonance at optical depth tv (defined by p(tv) = pv). 
(2) At the resonance, the X-mode and O-mode intensities are converted using eqs. (1331 and 1341 . (3) Integration of eq. @ is 
continued to Tmin- We use an analagous procedure for downward integration from Tmin to Tmax- 



3.1.5 Temperature correction procedure 

To integrate the radiative transfer equation, an initial temperature profile is assumed (the initial source function is set to 
Py/2). This initial profile is taken from a previously const ructed model with the same magnetic field and effective temperature, 
but without partial mode conversion (see|H2ALa?200?). In general, the solution to eq. using this profile will not satisfy 
eqs. To establish equilibrium, the initial temperature profile is corrected using the standard Unsold-Lucy procedure 

l)MihalaJIl973^ . The entire process is iterated until the deviations from radiative equilibrium, constant flux, and the relative 
size of the temperature correction are all less than a few percent. During a given iteration, the specific intensity calculated 
from the previous iteration is used to determine the source function. Thus, the source function must also converge to yield 
a self-consistent solution. Numerically, we find that the source function converges more rapidly than the other quantities 
consid ered above. For a m ore detailed discussion of the construction of self-consistent atmosphere models, see IHo fc Lail 
l)200ll) and lMihalaj l)l978^ . 



3.2 Partial mode conversion using photon Stokes parameters 

While the treatment described above captures the essential physics of the transfer problem, it is important to compare it to the 
exact solution obtained from integration of the transfer equations for the radiation Stokes parameters. As discussed in ij2l near 
the vacuum resonance, the modal transfer equation (|SJ breaks down because of the violation of the Faraday depolarization 
condition and collapse of the photon modes (see Figs. 4-5 of iLai fc Ha ii2003a) for the precise condition). 
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The radiation transfer equations for the Stokes parameters are given by JLai fc Hol l2003^: 
dl 



dr 



M I-S^ 
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(39) 



(40) 



(41) 



where cth = e' cos^ 6k b + r)' siv? 6kB — a, ai2 — igcosOks, 0-22 = e' — a — msin^OkB, Uafir — ?R.e{aai3), (Tapi — 9m((Ta/3), 
and /= (/ii, /22, f/f, K,) + , with /n = (Jj, + Q^)/2, /22 = (/i/ — Qv)/2- Note that eq. 1391 ignores scattering; the "em" sufBx 
on the source functions implies that terms proportional to 7re or 7r^ should be set to zero, as they are related to scattering 
contributions. The scattering contributions to eq. (I39II are derived in iLai fc Hj J2003ah . 

Away from the resonance, the modes discussed in H3.1.2I are well defined and are readily calculated from the Stokes 
parameters. Neglecting dissipative terms in the dielectric tensor, the transverse part of the mode polarization vectors can be 
written [see eq. 11711 1 



-i sin 6m, cos f 



(42) 



where 6m is the "mixing angle" defined by cos Um 
can be calculated from the Stokes parameters via 

= i [/. ± (cos 2emQ,. + sin 26»„K)] . 
Conversely, given the mode intensities, the Stokes parameters can be calculated using 



Kjf. I 1 + , sin 6'm = 1/ -y/ 1 + K']^ . The intensities of the ± modes 

(43) 



C + I^, (44) 

cos 26l„ (/+-/-) -2 sin 25™ (/+/-)i/^ cos A<^, (45) 

-2(/+/-)'/%inA0, (46) 

sin 25^(1+ - /") + 2 cos 2e„ (/+/-) '/"^ cos A,^, (47) 

where A0 — A(j)i + (oj/c) J^{n+ — n-)dz is the phase difference between the + and — modes. Note that Acj) is unknown, 
since the initial phase difference A^i between photons in the X and O-modes is random. To correctly evaluate the Stokes 
parameters from the specific mode intensities, one should sample A^ from a random distribution, and average over the results. 
Practically, we note that while the choice of A4) affects the values of the Stokes parameters, it does not change the specific 
mode intensities calculated from eq. (1431 . Therefore, the phase difference is unimportant for the comparison of the mode and 
Stokes parameter transfer equations (see H3.3II . 

In principle, eq. 1391 can be integrated from Tmax to Tmin using the initial condition I{Tmax) ~ (-8^/2, -8^/2, 0, 0)^ as in 
H3.ll However, this approach runs into a numerical difficulty: away from the vacuum resonance, differences in the indices of 
refraction for the two modes manifest as rapid oscillations in Q^, Uv, Vv, which are difficult to handle numerically. Thus, the 
direct solution of eq. 1391 over the entire range of integration is impractical. It is possible, however, to integrate eq. 1391 for 
a small range of t around the resonance. Using eqs. 1431 and 1441 - 1471 . we can quantitatively compare the result of such an 
integration with that obtained using the method of >I3.1.4I and thereby confirm the accuracy of the latter method (see ii3.31 . 



3.3 Numerical comparison between mode and Stokes transfer equations 

We consider a typical case, the propagation of a photon, initially polarized in the — mode, with energy E = 1.0 keV, 
propagation angle 6kB = 7i"/4, and magnetic field B = 10^* G. The temperature profile is held constant at T = 5 x 10^ K 
[eqs. I|44^ - I|47j are used to set the initial conditions for 139|l ]. Figure |21 shows the Stokes parameters as a function of optical 
depth near the resonance. These are obtained by integrating eq. I|39^ . The corresponding mode intensities are then calculated 
using eq. I|43|l and depicted in the top panel (solid and dashed fines). The dashed-dot and dotted lines show the results 
obtained from the integration of the mode equations with partial mode conversion [eqs. 1331 - 1341 ]. Note that the curves agree 
exactly except near the resonance where the modes are not well-defined. 

We have carried out many similar comparisons between the mode equations and the Stokes transfer equation. The close 
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Thompson optical depth t 



Figure 2. Evolution of the Stokes parameters across the vacuum resonance obtained by integrating the transfer equation 1391 . The 
parameters are B = 10^'^ G, Of^B = ■'r/4, T = 5 X 10^ K, and E = 1 keV. At high optical depth, the photon is in the — mode. In the top 
panel, the dotted line and the dot-dashed line show the mode intensities obtained using the method described in H3.ll [i.e.. solving the 
transfer equation j^J based on photon modes, but taking account of partial mode conversion through eqs. 03- J^], while the solid and 
dashed lines give the mode intensities based on the evolution of Stokes parameters (panels 2-4). Note the close agreement everywhere 
except near the resonance where the modes are not well-defined. 



agreement between the two methods establishes the validity of our method described in H3.ll i.e., integrating the mode eq. ||SJ 
and taking partial mode conversion into account using eqs. 13311 - 134(1 . 



3.4 Numerical grid and test of accuracy 

In solving eq. ®, we set up grids in Thomson optical depth, temperature, density, energy, and angle. The grid in optical 
depth {Td : d — 1, . . . , D} is equally spaced logarithmically with 15 — 20 points per decade (ppd). As discussed above, this 
grid spans 5 — 8 orders of magnitude to insure that photons are generated at densities higher than the X-mode decoupling 
depth, and that they are fully decoupled from the matter at the outermost layer. 

Care must be used in defining the energy grid. As mentioned before, vacuum polarization introduces a narrow spike in 
the X mode opacity. A prohibitively high e nergy grid resolu tion is required to properly resolve this feature. An alternative is 
to use the equal-grid method described bv lHo fc Lail i2003l) . In this case, each point of the energy grid is chosen to be the 
vacuum resonance energy [given by eq. correspo nding to a point on the optical depth grid: {En = Ev{t„) : n — 1, . . . , D}. 
This insures that the vacuum resonance is resolved. IHo fc Lail |20o3) point out that in the "no conversion" limit, this leads to 
an over-estimate of the integrated optical depth across the vacuum resonance. It is therefore important to investigate what 
effect this has on the emergent spectra. 

Figure |3| illustrates the effect of grid resolution on the spectra for two of the models presented in ^ The top panel 
shows the model with B = 10^* G, T^ft = lO'' K, which includes vacuum polarization in the opacities, but neglects the mode 
conversion effect (i.e., Pjump = 1). Over-estimation of the vacuum resonance in the X-mode opacity is expected to be strongest 
for this model, since modification of the emission spectrum is due solely to the enhaced opacity. The difference between models 
at 15 ppd and 20 ppd is negligible. Even at 10 ppd the difference is small, occuring mainly around the proton cyclotron line, 
as expected. 

The bottom panel of Fig. El shows the model with B = 10^* G, and Teg = 5 x 10® K, which includes partial mode 
conversion. At higher effective temperatures, the optical depth across the vacuum resonance becomes much greater than 
unity, thus, the error due to finite grid size becomes even less important. The lower grid resolution (10 ppd) model shows 
negligible deviation from the higher-grid resolution (20 ppd) model, even around the proton cyclotron feature. This behavior 
is typical of all models with _B = 5 x 10^^ G. 
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Figure 3. Spectra showing the effect of grid resolution on thermal emission. The top panel shows hydrogen atmosphere models with 
B = lO^-* G, Tcff = 10^ K, which include vacuum polarization but neglect mode conversion. As the number of grid points per decade 
(ppd) are increased, the curves quickly converge to the 20 ppd case. At low resolution the error mainly occurs around the proton 
cyclotron feature, and is negligible elsewhere. The bottom panel shows hydrogen models with B = 10^* G, Tgff = 5 X 10^ K, which 
include vacuum polarization and mode conversion. At higher effective temperature, the difference between models with varying grid 
resolution is negligible. 



4 RESULTS: ATMOSPHERE STRUCTURE, SPECTRA, AND EMISSION BEAM PATTERN 

We now present the results of our atmosphere models. We consider B = 4x 10^^, 7x 10^^, lO^"*, 5x 10^" G, and Tofi = lO*", 5x 10' 
K, for both H and He compositions. 



4.1 Atmosphere structure 

Figure0|shows the temperature profile as a function of Thompson optical depth r for the H atmosphere model with B — 10^"* 
G, and Tcff = 10® K. To understand the effect of vacuum polarization, we show the results based on four different ways of 
treating vacuum polarization; (1) vacuum polarization effect is completely turned off ("no vaccum"); (2) vacuum polarization 
is included, but the mode conversion is neglected (Pjump ~ 1, "no conversion"); (3) vacuum polarization is included, and 
complete mode conversion is assumed (Pjump = 0, "complete conversion"); (4) vacuum polarization is included with the 
correct treatment of the resonance ("partial conversion"), using Pjump calculated from eq. I^J- We see that models which 
include vacuum polarization show higher temperatures over a wide range of r for the same Tcff than models which ignore 
vacuum effects. This temperature increase is due to the X-mode opacity feature at the vacuum resonance (see Fig. 0. In 
general, atmosphere structure is determined by the radiative equilibrium condition, and inspection of the individual terms of 
eq. ||SJ reveals how the resonance affects the temperature profile. The mode absorption opacities obey the relation no 3> kx 
except at the energies E — EBi,Ev- Thus, kx{Bi,/2) can be neglected relative to k.o{Bi,/2) in eq. (|HJ. In the absence of 
vacuum polarization, the O mode largely determines the atmosphere structure due to the weak interaction of X mode photons 
with the medium. However, when the resonance spike in the X mode opacity is present, i^xJ^ cannot be neglected relative to 
hqJ!^', in fact, this occurs over a large bandwidth for which ^ . The result of this enhanced interaction is to increase 
the overall temperature. Adding the effect of mode conversion further increases the temperature over a large range of optical 
depth. This is due to heat deposited by converted X mode photons, which interact with the large O mode opacity after passing 
through the vacuum resonance. The temperature profile for the partial mode conversion model (shown by the solid curve of 
Fig. |1J closely follows the result for the no conversion model (shown by the dotted curve) for the small optical depths at which 
low energy photons decouple. This is because for these photons, E < Ead, is satisfied and mode conversion is ineffective. For 
larger optical depths, at which higher energy photons decouple, E > Ead, and mode conversion is more effective, thus the 
"partial conversion" result lies between the "no conversion" and "complete conversion" limits. 

Figure 13 shows the temperature profile for the B = lO'^* G, Teff = 5 x 10® K model. This higher temperature model 
shows the same basic features as the low- Teff model in Fig. 2] In this case, the energy flux is carried by photons with higher 
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Figure 4. Temperature profile for the hydrogen atmosphere model with B = 10^"^ G, T^ff = 10® K. The four curves correspond to 
different ways of treating the vacuum polarization effect: (1) no vacuum (short-dashed curve); (2) no conversion (dotted curve); (3) 
partial conversion (solid curve); and (4) complete conversion (long-dashed curve). 




Figure 5. Same as Fig. El except for Teff = 5 X lO'' K. 



energies, and the adiabatic condition [E > Ead) is more readily satisfied, leading to effective mixing of photon modes. Thus 
we see that at large optical depths (r > 0.1), the partial conversion profile closely follows the complete conversion curve. At 
lower optical depth, the partial conversion profile lies between the complete conversion and no conversion curves. 

Finding self-consistent temperature profiles is the most time consuming step in atmosphere modeling. Once the profile is 
known, the emergent radiation can be obtained by a single integration of the transfer equation. To facilitate future work on 
NS atmospheres and related applications, we provide fitting formulae for the models presented in this paper. Formulae are 
provided only for models incorporating vacuum polarization with partial mode conversion. The fits are vafid over the optical 
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Table 1. The parameters in the fitting formulae [Eq. 1481 1 for the temperature profiles for atmosphere models with different magnetic 
field strengths, effective temperatures and compositions (ionized H or He) 
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depth range r = 10 — 2 x 10 . Each model is fit by the function 



FT r ^1 - J ai + a2 Ax + as Ax'-^ + a4 Ax^ + as Ax'^ + ae Ax^ Tmid < r < 2 x 10"*, 
log,o[U(r)\ - i a^ + a2Ax + hAx^ + b4Ax-' + hAx^ + hAx^ IQ-' < r < t^,^, 



where x = logj^Q(r), A and Tmid denotes the break between the two parts of the fit necessary to describe the 

temperature profile. The parameters for each model are summarized in Tabled 



4.2 Spectra 

Figure |S| presents the spectrum for the hydrogen atmosphere model with B = 10^"* G and Toff = 10^ K. The results for the 
four different ways of treating vacuum polarization (see H4.1II are shown. These spectra correspond to the temperature profiles 
depicted in Fig.|l] 

When the vacuum polarization effect is neglected, the spectrum of a magnetic, ionized H atmosphere model generally 
exhibits two characteristics: (1) a hard spectral tail (compared to blackbod y) due to the non-grey free- free opacity (k^ decreases 
with increasing photon energy; IShibanov erai]ll99i IPavlov et alJll995^: (2) a significant proton cy clotron absorption line 
when Ebz is not too far away from the blackbody peak (~ SknT) l|Zane et alJboOll: IHo fc LailEoOlll . Vacuum polarization 
tends t o soften th e hard spect ral tail and suppr ess (reduce) the proton cyclotron line. These effects are discussed extensively 
in §4 of lHo fc Lai (|200,'^ ) . and lLai fc Hoi (|22q3)- In Fig- El of the models that include vacuum polarization effects display 
a large reduction in the equivalent width (EW) of the proton cyclotron feature due to the modification of the temperature 
profile by the vacuum resonance (see H4.1II and the mode conversion effect. The spectra also show softening of the hard spectral 
tail relative to the no vacuum case, though they are all still harder than blackbody. The "partial conversion" curve appears 
as an intermediate case between the "complete conversion" and "no conversion" extremes. The adiabatic regime where mode 
conversion is efficient is clearly visible: for E > Ead ~ 2 keV, the "partial conversion" curve begins to follow the "complete 
conversion" curve. 

This transition from "no conversion" to "complete conversion" is further illustrated by Fig. |7| which shows, for several 
photon energies E (and a given direction of propagation, 9kB) the evolution of the specific mode intensities as a function of 
optical depth. At the resonance depth (denoted by the vertical lines), the X mode photons encounter the vacuum induced 
spike in opacity, and mode conversion occurs [governed by the non-adiabatic jump probability Pjump; see eq. Q]. As the 
energy is increased (from the second to bottom panels), mode conversion becomes increasingly more effective. 

Figure I2I shows the spectrum for the B = 10^* G and Tcft = 5 x 10^ K hydrogen model. All the calculations that include 
vacuum effects show strong suppression of the ion cyclotron feature and significant softening of the hard spectral tail. At higher 
effective temperatures, there is a smaller difference between the no conversion, partial conversion, and complete conversion 
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E (keV) 

Figure 6. Spectra for hydrogen atmosphere models with B = 10^* G, Teff = 10® K. The four cases described in the text are shown: 
(1) no vacuum (short-dashed cuve); (2) no conversion (dotted curve); (3) partial conversion (solid curve); and (4) complete conversion 
(long-dashed curve). The light dashed-dot curve shows the blackbody spectrum with T = 10^ K. For all three cases that include vacuum 
effects, the proton cyclotron feature is strongly suppressed and the high energy tail is softened relative to the no vacuum case. The 
"partial conversion" curve is seen to be intermediate between the "no conversion" and "complete conversion" limits. 




Thomson optical depth t 



Figure 7. The transition from "no conversion" to "complete conversion" for a photon propagating in a hydrogen atmosphere with 
B = 10^* G, Tgfj = 10® K. The top panel shows the temperature profile for this model. The bottom three panels show the evolution of 
the specific mode intensities for energies E = 0.35, 0.65, 2.72 keV and B^g = it/ 4. The dashed line shows the X-mode intensity, the solid 
line the O-mode intensity, and the dotted vertical lines specify the location of the vacuum resonance. For E = 0.35 keV: Hp = 0.69 cm, 
E < -Eadi and Pjump = 1-0, leading to minimal mode conversion. For E = 0.65 keV: Hp = 0.78 cm, E ~ Eadi ^'Hd Pjump = 0.65, leading 
to partial mixing of the modes. For E = 2.72 keV: Hp = 0.98 cm, E > Ead, and Pjump = 0.15, leading to nearly complete conversion of 
the modes. 
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E (keV) 

Figure 8. Same as Fig.|H] except T^ff = 5 X 10^ K. The strength of the proton cyclotron feature is strongly suppressed for models which 
include vacuum polarization. At E > 5 keV, the models which include mode conversion are softer than those without, though all the 
atmosphere models are still harder than blackbody. 
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Figure 9. Spectrum of the hydrogen atmosphere model with B = 4 X 10^^^ G, TcS = 10^ K, calculated for two cases: partial mode 
conversion (solid curve), and no vacuum (dotted curve). The light dashed line shows the blackbody spectrum with T = 10® K. 



cases. In this regime, the optical depth across the vacuum resonance is much greater than unity, and the decoupling o f X-mode 
photons occurs at the resonance density whether or not mode conversion is taken into account (see lLai fc HolEoOSi) . At high 
energies {E > 5 keV), the models which include mode conversion are softer than those which do not. 

Figures f9ll2l deDict atmosphere models at Tcff = 10® K with varying magnetic field strengths, comparing the "no vacuum" 
and correct "partial conversion" results. For B = 4 x lO" G < Bi (Fig.EJ, the vacuum resonance lies at a lower density than 
the X-mode and 0-mode photospheres, and vacuum polarization has a negligible effect on the spectrum, reflected by the close 
agreement between the "no vacuum" and "partial conversion" curves. For the hydrogen atmosphere model with B = 7 x 10^^ 
G (Fig, nut , and the helium atmosphere model with B = 10^* G (Fig. 111^ . vacuum polarization affects the spectrum. For 
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E (keV) 

Figure 11. Same as Fig. |2l except for lielium composition with B = 10^* G. The ion cyclotron line width is reduced by vacuum 
polarization, though it has a larger equivalent width than the model for hydrogen, due to the location of the line near the maximum of 
the continuum emission. 



these models, B > Bi, and the ion cyclotron line lies near the blackbody peak. Thus, it is particularly important to treat the 
vacuum resonance correctly, taking partial mode conversion into account to calculate the line width. For the B — 5 x lO^'' 
G model (Fig. the spectral feature at Esi is outside the energy band of observational interest. We note that at such a 
high fie ld and low eff ective temperat ure, the atmosphere should contain a significant fraction of bound atoms and molecules 
JHo et al . 2003 ; Potekhin et alj|2004l) , so the fully ionized model shown in Fig. I12l is not realistic. 

Figures ll 31141 show atmosphere models at magnetic field strength B = 5 x 10^^ G and Tcff = 5 x lO" K, for hydrogen and 
helium compositions. At this effective temperature, the ion cyclotron feature lies close to the blackbody peak, and the effects 
of vacuum polarization on the line width and the spectral tail are rather pronounced. 
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' B = 5x101" G 




E (keV) 

Figure 12. Same as Fig. CHI except with B = 5 X IQl"* G. 



B = 5x1014 G 




E (keV) 

Figure 13. Same as Fig. 1121 except for T^g = 5 X lO'' K. The proton cyclotron is strongly suppressed by vacuum polarization. The hard 
spectral tail is softened considerably relative to the no vacuum case, though it is still harder than blackbody. 

4.3 Emission beam pattern 

Calculations of observed NS lightcurves and polarization signals (see ij^J are critical for interpreting observations. An important 
ingredient of such calculations is the angular beam pattern of surface emission. Figure [T^ shows the radiation intensity from 
a local patch of NS atmosphere (for the Teg = 10*' K models presented in >I4.21 as a function of emission angle relative to the 
surface normal at several photon energies. The heavy curves show models that include vacuum polarization effects, while the 
lightcurves show models that neglect vacuum polarization. Magnetized atmosphere models which neglect vacuum polarization 
have a distinctive beaming pattern, consisting of a thin "pencil" feat ure at low emission angles and a broad "fan" beam at 
large emission angles, with a prominent gap between them (e.g., c.f.. lOzelboOlh . This gap tends to increase with increasing 
photon energy. The detailed shape of the emission beam pattern is determined by the temperature profile and the anisotropy 
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B = 5x1014 G, He^.' 




E (keV) 

Figure 14. Same as Fig. 1131 except for helium composition. Vacuum polarization strongly suppresses the ion cyclotron feature, and 
softens the hard spectral tail relative to the no vacuum case. 




20 40 60 80 20 40 60 80 
^'em (degrees) 

Figure 15. The specific intensity from a local patch of NS as function of emission angle 0em (= cos~^ /i, the angle between k and the 
surface normal) for hydrogen models with T^g = 10^ K, and B = 4 X 10^^ G, 7 X lO^^ G, IQl* G, 5 X IQl* G. The heavy curves show 
models that include vacuum polarization effects, while the lightcurves show models that neglect vacuum polarization. The latter models 
exhibit a characteristic beaming pattern, with a thin "pencil" shape at low emission angles and a broad "fan" at large emission angles. 
Inclusion of vacuum polarization tends to reduce the gap and lead to a broad, featureless beaming pattern at large magnetic fields. 



of the mode opacities. As shown by Fig. 1151 vacuum polarization tends to smooth out the gap, leading to a broad, featureless 
beam pattern at large magnetic fields. The broadening of the beaming pattern is due to the alteration of the temperature 
profile by the spike in opacity at the vacuum resonance and the mode conversion effect. 

Figure ITfil shows the specific radiation intensity from a patch of NS atmosphere (for the Tofi — 10^ K models presented 
in H4.2t at several emission angles. We see that the shape of the spectrum and EW of the absorption feature can change 
significantly depending on the emission angle and whether or not vacuum polarization effects are included in the calculation. 
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E (keV) 

Figure 16. The specific radiation intensity from a local patch of NS atmosphere at several emission angles for hydrogen models with 
Tcft = lO'^ K, and B = 4 X 10^^ G, 7 X lO^^ G, 10^* G, and 5 x 10^* G. The shape of the spectrum and EW of the absorption feature 
change significantly with emission angle. 



5 POLARIZATION OF THE ATMOSPHERE EMISSION 

Thermal emission from a magnetized NS atmosphere is highly polarized. This arises from that fact that the typical X-mode 
photon opacity is much smaller than the O-mode opacity,^ kx ~ {E /Ebc)^ no ^ hq. Thus, X-mode photons escape from 



hieh degree fas hieh as 100%) (see. e.e..lGnedin & Sunvaevlll974: iMeszaros et alJIlQSa: 


IPavlov & ZavhnlboOOl) 


There has been some recent interest in X-ray polarimetry for NSs 


ICosta et alj 


|200l|, 


boodl. Observations of X 



polarization, particularly when phase-resolved and measured in different energy bands, can provide useful constraints on 
the magnetic field strength and geometry, the NS rotation rate and compactness. This will be highly complimentary to 
the information obtained from the spectra and lightcurves. Moreover, as we show below (see also Lai & Ho 2003b), vacuum 
resonance gives rise to a unique polarization signature in the surface emission, even for NSs with "normal' {B ~ 10^^ - lO" G) 
magnetic fields. 

Below, we calculate the observed X-ray polarizaztion signals in the case when the emission comes from a rotating 
magnetic hotspot on the NS surface. Although this represents the simplest situation, it captures the essential properties of 
the polarization signals, and the result can be carried over to more general situations. See the end of >l5.3l for a discussion of 
the case when emission comes from an extended area (or the whole) of the stellar surface. 



5.1 Geometry and lightcurves 

To calculate the observed lightcurves and polarization signals, we set up a fixed coordinate system XYZ with the Z-axis 
along the line-of-sight (pointing from the NS toward the observer), and the X-axis in the plane spanned by the Z-axis and 
fl (the spin angular velocity vector). The angle between ft and ez is denoted by 7. The hotspot is assumed to be at the 
intersection of the NS surface and dipole magnetic axis, which is inclined at an angle rj relative to the spin axis. As the star 
rotates, the angle O between the magnetic axis fj, and the line of sight varies according to 

cos 6 — cos 7 cos rj — sin 7 sin rj cos 1/) , (49) 



° This applies under typical conditions, when the photon energy E is much less than the electron cyclotron energy Eg^ , is not too close 
to the ion cyclotron energy Egi, the plasma density is not too close to the vacuum resonance (see the text), and 9f^g (the angle between 
k and B) is not close to 0° or 180° . 
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Figure 17. Lightcurve and polarization as a function of rotation phase for a NS hotspot with B = 10^^ G, T^ff = 5 X 10^ K. The angle 
of the spin axis relative to the line of sight is 7 = 30°, and the inclination of the magnetic axis relative to the spin axis is ?7 = 70°. Note 
that the sign of the Fq Stokes parameter is opposite for low and high energy photons; this implies that the pianos of polarization for 
low and high energy photons are perpendicular. This is a unique signature of vacuum polarization for models with B < Bi. 



w here il} = (Ot + consta nt) is the rotation phase {tp — when /x lies in the XZ plane). We use the simplified formalism derived 
by iBeloborodovl i2002l) to calculate the observed spectral fiux from the area dS of the hotspot, which takes the form: 

Fobs ^i}-^) ^4^cm)cos6i„„-^, (50) 

where r-g = 2GM/c? is the Schwarzschild radius, R is the NS radius, and ^cm (the angle between the photon direction and 
the surface normal at the emission point; /i = cos^em) is related to through: 

COSee^ = ^+(l-^)cOSe. (51) 

For R > 3rg, the spectral fiux calculated using eq. 15U1 differs from the exact treatment fsee IPechenick et al.lll98^ by less 
than 1%. 

The top panels of Figs. ll7ll9l show lightcurves for NS models with several magnetic field strengths, at a range of energies 
= 0.5 - 3 keV, with geometry 7 = 30°, r? = 70°. 



5.2 Observed linear polarization signals 

The atmosphere models presented in JSl^nd ?S]give the specific intensities /^(^om) and /°(Scm) of the two photon modes, 
emerging from the NS atmosphere (outside the vacuum resonance layer). To determine the observed polarization signals, it is 
important to consider propagation of the polarized radiation in the NS magnetosphere. In the X-ray band, the magnetospheric 
dielectric properties are dominated by vacuum polarization (Heyl & Shaviv 2002). Heyl et al. (2003) (see also Heyl & Shaviv 
2002) evolved the Stokes parameters along photon geodesies in the magnetosphere and showed that the observed polarization 
is determined at the so-called "polarization limiting radius," at a distance far from the NS surface. Below we present a simple 
calculation of the propagation effect and observed linear polarization (see also Lai & Ho 2003b) . 

Consider a photon emitted at time ti from the hotspot, with rotation phase xpi — ilti. The emission point has polar angle 
Oi (relative to the fixed XYZ frame) given by eq. 14911 with = il^i, and azimuthal angle (pi given by 

sin ri sin ipi 

tanifi — ; . (52) 

sm ri cos 7 cos Tpi + cos ri sm 7 

This is also the angle [ipi = lpb{R)] between the projection of the magnetic axis in the XY plane and the X-axis. After the 
photon leaves the star, it travels towards the observer, with a trajectory given by 

r = (i?sineiCos</)i + Ax) X + (i?sin 6; sin 0i + Ay)Y + iRcosQ^ + s + Az) Z, (53) 
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Figure 18. Same as Fig. 1171 with B = 7 X 10^'' G. For this model, B ~ Bi, representing the transition between the models shown in 
Figs, nn and [T^ 




Figure 19. Same as Fig. 1171 with B = 5 X 10^'* G. At this field strength, B > B;, and the vacuum resonance lies between the O and 
X-mode photospheres. Thus, the sign of the Fq Stokes parameter for low and high energy photons is the same. Note that in the top 
panel, the flux for the 0.5 keV case is multiplied by a factor of 10 relative to the other curves. 



where s — cAt = c{t — ti), Ax,y,z are relativistic corrections (which, as we will see shortly, are unimportant for the polarization 
signals), and X, Y, Z are unit vectors. As the photon propagates in the magnetosphere, it will "see" a changing stellar 
magnetic field, given by B = — V(/x ■ r/r^), where'^ 

(J, ~ ^ [(sin rj cos 7 cos -ip + cos rj sin 7) X + sin rj sin rpY + (cos r] cos 7 — sin rj sin 7 cos 1/;) ^] , (54) 
with tp — Qt = tpi + flAt — ipi + s/ri (here ri = c/Q is the radius of the hght cyclinder). The photon's polarization state 

^ We restrict the propagation to the "near zone" of the star, i.e., r < r; = c/Q. 
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will evolve adiabatically following the varying magnetic field that the photon experiences, up to the polarization limiting 
radius rpi, beyond which the polarization is frozen. Since we anticipate rpi 3> R, we only need to consider the region far 
away from the NS: for r ^ R, the photon trajectory is simply r ~ s Z, and the magnetic field along the photon path is 
B ~ {2fj,zZ — fixX — /iyy)/r"^, with r ~ s. This magnetic field has a magnitude 

B(s) = (^ — ^ + 3 (cost; cos 7 — sin 7 sin?) cos i/j)^] ^''^ , (55) 

where Bg — 2fj,/R^ is the magnitude of the (dipole) surface field at the magnetic pole. The magnetic field is inclined at an 
angle 9kB to the line of sight, and makes an azimuthal angle ifiB in the XY plane such that: 

.2^ , . 1 - (cos r; cos 7 - sin 7 sin 77 COS ■0)^ 

sm OkBis) = — — : : -2, (56) 

1 + 3 (cos rj cos 7 — sm 7 sm ri cos tp) 

I \ sin 77 sin V' 

tan V3s(s) = : —. -. (57) 

cos r\ sm 7 + cos 7 sm r\ cos '^) 

Recall that in Eqs. 1551 - 157j . tp — tpi + s/ri — Qti + s/vi and s is the affine parameter along the ray. 
The wave equation for photon propagation in a magnetized medium takes the form 
2 

V X (/i ■ V X E) = ■ E, (58) 

where E is the electric field (not to be confused with the photon energy), and e, fi are the dielectric and inverse permeability 
tensors, respectively. In the magnetized vacuum of the NS magnetosphere, they are given by e = al + qBB and p, — al + mBB. 
Solving eq. ll^ for EM waves with E oc e'*"'"''^* yields the two modes (in the XY basis) 

eo = (cos(^B,sinv7B), ex = (— sin (^s, cos ^s), (59) 

with indices of refraction no — l + {q/2) sin^ dks and nx — 1 — ('ti/2) sin'^ dks ■ A general (transverse) EM wave can be written 
as a superposition of the two modes: 

E = ^oeo+^xex. (60) 

Following the steps of Lai & Ho (2002) [see their eq. (15)], we derive the following equations for the evolution of the mode 
amplitudes: 



1 



■A'o \ ^ ( -(^/c)An/2 i^'s \ ( -^o 
A'^ J I -iip's {Lo/c)An/2 1 \ Ax 



(61) 



where the prime (') denotes a derivative with respect to s, and An — no — nx = ^{q + m) sin^ 9kB- The condition for adiabatic 
evolution of photon modes is 

((.j/c)An> 2v9s. (62) 

Near the star (r ~ R), ips ~ 1/r, and the adiabatic condition is easily satisfied. Far from the star (r ^ R), using eq. 1551 . we 
write: 

^"=3^(:fe) -n^^.B=9.94xlO-«S?2(^)%s, (63) 
where 

Fb = 1 — (cos r; cos 7 — sin 7 sin 7; cos , (64) 
and B12 = Bs/(W^^G). From eq. we have 

—T— = —F^, (65) 
ds ri 

with 

F^ = (sin'^ rj cos 7 + sin 77 cos rj sin 7 cos t/") / Fb ■ (66) 
The polarization-limiting radius r-pi is set by the condition uj^njc = 2Lpg. Substituting in eqs. 16311 and 1651 we find * 

* Our expression for rpi differs from that given in Heyl &; Shaviv (2002) and Heyl et al. (2003). 
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where /i is the spin frequency n/(27r) in Hz, and Fb, are slowly varying functions of phase and are of order unity. Note that 
the above analysis is valid only if rpi ti/2, since beyond the light-cylinder radius the magnetic field is no longer described 
by a static dipole. Thus we require that 

^ ^ 6.84 X 10-^ ( El^fA R,, ^ 0.5, (68) 
ri \ F^ J 

where Rio = R/{10 km). 

Beyond rpi, the polarization state of the photon is "frozen." Using eq. 1591 1. the observed Stokes parameters in the observer 
coordinate system (XYZ) are given by 

/ = /o + /x, (69) 
Q ~ (Jo -/x)cos2v9B(rpi), (70) 
U ~ (Jo-/x)sin2<^s(rpi), (71) 

where Jo oc /°(Scm) and Jx oc /^(^cm) are the specific mode intensities emitted at the NS surface [calculated with our 
models described in 33]and ^ and corrected for the general relativistic effect; see eqs. 15UII - 15HI ]. and v'-b(''pi) is evaluated at 
s ~ r = Tpi. From eqs. 11521 and 1571 . with '(/'(j'pi) = tpi + rpi/ri, we see that the effect of NS rotation is to shift the polarization 
lightcurve by a phase rpi/rj. For slow rotation r^i/ri <^ 1 [see eq. 1681 1. this shift is small, and ^PBifpi) ~ ipi + n. We calculate 
the observed spectral fluxes Fj — F, Fq, Fu associated with the intensities /, Q, U using the standard procedure described 
in 

The middle and bottom panels of Figs. ll7ll9l show the phase evolution of the Stokes parameter Fq and the degree of 
linear polarization, both normalized to the observed spectral flux. Note that Q is defined such that Q = 1 corresponds to 
linear polarization in the plane spanned by the line of sight Z and the NS spin axis. In Fig. 1171 we consider emission from a NS 
hotspot with B — 10^^ G and T^g = 5 x 10^ K. Note that the value of Fq for low energy photons {E ^ 1 keV) is of opposite 
sign to that of high energy photons (i? ^ 3 keV). This implies that the planes of polariza tion for low and h igh energy photons 
are perpendicular. This is a unique signature of vacuum polarization first identified bv lLai fc Hoi ll2003lJ) . which occurs for 
B < Bi [see eq. because the vacuum resonance appears outside the O-mode photosphere. Below the vacuum resonace layer 
(p > pv), the X-mode flux dominates over the O-mode. For low energy photons, E ^ Ead, mode conversion is inefficient, and 
the emergent flux is dominated by the X-mode; for high energy photons E ^ Ead, mode conversion is efficient, rotating the 
plane of polarization, and the emergent flux is dominated by the O-mode [see Fig. 2 of Lai & Ho (2003b)] . 

Fig. 1191 shows the same result for the model with B = 5 x 10^'' G, Tcs — 5 x 10^ K. In this case, B > Bi, the vacuum 
resonance appears inside the O and X-mode photospheres, and the emergent radiation is always dominated by the X-mode. 
As expected, the planes of polarization for low and high energy photons are aligned. Fig. 1181 shows an intermediate case, 
where B ~ Bi — to calculate the spectra and polarization signals for such models, it is particularly important to incorporate 
partial mode conversion properly. The distinct behavior between the low field and high field cases is illustrated by Fig. 1201 
which shows the phase-averaged Fq Stokes parameter as a function of photon energy for several values of the magnetic field 
strength. The low-field cases show the characteristic rotation of the plane of polarization between low-_E and high-_E, whereas 
the high-field cases do not. 

Note that in the above analysis, the observed linear polarization fraction Hl (the bottom panel of Figs. 1171191 is the 
same as the value just outside the emission region, |Hem|, i.e.. 

The polarized fluxes are simply 

Fq ~ F/Hem cos2((9fl(rpi), Fc/ ~ F/Hem sin 2(^3 (rpi). (73) 

In i|5.3l we shall see that for rapidly rotating NSs, the observed Hl will be somewhat smaller than jHemj because of the 
generation of circular polarization around rpi. 



5.3 Circular polarization 

Circular polarization of surface emission may be generated for sufficiently rapid rotating NSs due to the gradual wave mode 
coupling/decoupling around rpi. While the linear polarization signals can be adequately described and calculated using the 
simple analysis given in 35.21 to calculate the circular polarization, quantitative solutions of the evolution equations for the 
modes or the Stokes parameters in the magnetosphere are necessary. 

Consider the mode evolution equations 1611 . For given initial values (e.g., mode amplitudes at r; <^ rpi), the solution of 
the equations depend on two parameters: C, defined by (Lj/c)An « C/r^ and ip'g (since ip varies by a rather small amount 
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Figure 20. Phase average of Stokes parameter Fq as a function of plioton energy for a rotating NS hotspot at magnetic field strengths 
B = 10^^, 4 X 10^^, 7 X 10^^ and 5 x 10^** G, with T^g = 5 x 10^ K. Note that the sign of (Fq) changes between low and high photon 
energies for the low-fiold cases, corresponding to rotation of the plane of polarization. 

along the ray path, Fb and Fy, are nearly constant). Alternatively, since Tpi is determined by C/r® = 2(p'g, the solution 
depends only on the dimensionless parameter F, defined by 

r = rpi^'s = ^F^~ 6.84 x lO"" (SiBJzFsF^) 7?io (74) 



Indeed, if we define x — r/vpi, eq. 16H can be rewritten as 




Note that while eq. 16111 is valid for all radii, the magnetostatic approximation for the field of a rotating dipole breaks down 
beyond the light-cylinder radius. Since ip'g ~ !/»■;, eq. I|75|l is valid only for F ~ rpi/ri ^ 0.5. 
The Stokes parameters can be written in terms of the mode amplitudes as 



/ = \Ex\^ + \Ey\' = \Ao\' + \Ax\^ (76) 

Q = \Ex{^ -\Ey\^ =cos2<pb(\Ao\^ -\Axf)-2sm2,pB?fie{AoAi,), (77) 

U = 2?lie{ExEY)^sm2ipB{\Aof -\Ax\') + 2cos2ipB?Sie{AoA*x), (78) 

V = 2Qm{ExEY) ^2'Qm{AoAx) ■ (79) 



Figure 1211 gives some examples of the results of numerical integration of eq. H75^ . We start the integration at radius 
Xi = Ti/rpi such that the adiabatic condition is well satisfied (we typically choose Xi ^ 1/3). In these examples, the initial 
values are = 1, Ax = 0. After obtaining Ao{x) and Ax{x), we calculate the Stokes parameters using eqs. 17611 - 179(1 with 
(pB ~ 'f'si''' ~ ''0 = r(a; — Xi) (adding a constant to </5s will affect Q and U , but not V). We see that for x 1/2, the photon 
modes evolve adiabatically, and thus Q oc cos2</5s, U oc sin2<^s [see eqs. 170II - 17H ] and V ~ 0. Around x = 1 (r = rpi), 
the modes couple and circular polarization is generated. For x 2, the values of the Stokes parameters are "frozen" and no 
longer evolve. NSs with large F (corresponding to rapid rotations; e.g., / = 40 Hz, B — 10^^ G, -B = 1 keV yields F ~ 0.3) 
can generate appreciable circular polarization, with V{r <x>)/I ~ —14%. As the NS spin frequency decreases, the resulting 
\V/I\ decreases. 

Alternatively, using the relations 1761 - 179j and the mode evolution equation 161L we can derive evolution equations for 



the Stokes parameters: 

/' = 0, (80) 

Q' = (ojAn/c) sin2(/PB, (81) 

U' = -(a;An/c) cos2vPB, (82) 



© 2006 RAS, MNRAS OOO.ITllsTI 



Atmosphere Models of Magnetized Neutron Stars 25 




X 



Figure 21. Evolution of tlie radiation mode amplitudes (top panel) and Stokes parameters (bottom three panels). The solid lines are 
for r = 0.3 and the dashed lines for F = 0.1. The polarization limiting radius is shown as the vertical dotted line. The initial values (at 
a small x = Xi) are = 1> -Ax = 0, Q = / = l, (7 = and V = 0. At distances x ^ 0.5, the modes evolve adiabatically. At r ~ Tpi 
the modes begin to couple, generating appreciable circular polarization. At x 2, the values of the Stokes parameters are "frozen" . 



V' = ~{LuAn/c)Q sm2ifiB + {u}An/c)U cos2ipB. (83) 

Since the vacuum contribution to the dielectric tensor includes no dissipation, /' = as expected. Equations 18UII - 183(1 can be 
evolved numerically to calculate the observed Stokes parameters. Again we start the integration at a radius < rpi such that 
the adiabatic condition is well satisfied. Since the circular polarization does not depend on the orientation of the XY axes, 
we set the initial conditions at by rotating the coordinate system azimuthally so that 7 = 1, (7 = (this also corresponds 
to choosing = or (/3_b = at r = r;), and Qi = Hem ^ 1, where Hem is the linear polarization fraction just outside the 
atmosphere [see eq. 17211 1. Since the radiation emerging from the NS surface is linearly polarized, we set V{ri) = 0. Equations 
18UI I- 183|I are then integrated to a distance beyond rpi. 

For a given initial value of the linear polarization fraction Hem at a small <IC rpi, the solution of eqs. 181ll - 183|l depends 
only on the dimensionless parameter F = r^np'g. Again, if we define x = r/r^i, eqs. 181ll - 183|l can be rewritten as 

2V 

dQ/dx = — Vsin2^fl, (84) 

2r 

dU/dx = rrV cos2ipB, (85) 

X" 

2r 

dV/dx = (Q sin2(/Pi3 - 17 cos2(/Pi3) , (86) 

x" 

with ifB ~ F{x — Xi). We are interested in the value of V at x = Xf ^ 1. Figure l^ shows some examples of the integration of 
Eqs. 1841 - 186(1 . Not surprisingly, the results are in exact agreement with those obtained using the mode evolution equations. 

We have calculated the circular polarizations produced by rotating NSs with different values of F. Our numerical results 
[see Fig. I22ll 1 show that the generated circular polarization is given by the expression 

V/I ^ -0.60 He™ sign(^'s) Irpi ipB^^"^ = -1.5 x 10"=^ Uem (El Bf^ FbY''' fi F^. (87) 
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Figure 22. The magnitude of the observed circular polarization fraction \V\/I a,s a, function of F. The linear polarization fraction {Hem) 
just outside the atmosphere is assumed to be 100%. The dashed line depicts the fitting formula, eq. I87i . which agrees with the numerical 
solution to within 1% for F < 0.4. Note that solutions with F ^ 0.5 are incorrect since the magnetic field around the polarization limiting 
radius is no longer described by the near-zone field of a rotating dipole as adopted in our calculations. 



This expression is accurate to within one percent in the regime F ^ 0.4 [see Fig. 122ll 1. Recall that Fb ~ ~ 1, so eq. 1871 
provides a quick estimate for the magnitude of the circular polarization of NS surface emission. 

The above analysis shows that substantial circular polarization can be generated only for NSs with sufficiently rapid 
rotation and magnetic field strength. Given the photon energy, NS spin frequency, and dipole magnetic field strength, the 
degree of circular polarization from a NS can be calculated using eq. H87|l . 

Fig. 1231 shows the phase- resolved observed radiation Stokes parameters for a NS with B = 10^^ G and Tes — 5 x 10^ K, 
rotating at / = 50 Hz, with magnetic field and spin geometry 7 = 30°, rj = 70° (this is the case shown in Fig. I17II . The 
solid curves are numerical solutions to the Stokes parameter equations of transfer in the NS magnetosphere, while the dotted 
curves are calculated using the method described in § 15.21 The latter method assumes Fv ~ 0, but yields results for Fq and 
Fu that are quite close to those of the numerical integrations. For a rapidly rotating NS, substantial circular polarization is 
generated, with \Fv/Fi\ reaching 0.2 — 0.3 in the hard X-ray band (see also Fig. I22|l . 

Finally, we note that although the specific results presented in this section refer to emission from a hot polar cap on 
the NS, we expect many of our key results (e.g., rotation of the planes of linear polarization between E ^ 1 keV and at 
-E <; 4 keV due to vacuum polarization for B <, 7 x 10^^ G) to be valid in more complicated models (when several hotspots or 
the whole stellar surface contribute to the X-ray emission). This is because the polarization- limiting radius (due to vacuum 
polarization in the magnetosphere) lies far away from the star [see eq. 16 7| where rays originating from different patches of 
the NS experience the same dipole field (Heyl et al. 2003). Our results therefore demonstrate the unique potential of X-ray 
polarimetry in probing the physics under extreme conditions (strong gravity and magnetic fields) and the nature of various 
forms of NS. 
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Figure 23. Observed radiation Stokes parameters for a NS with B = 10^^ G, Teff = 5 X 10'^ K, / = 50 Hz, 7 = 30°, and rj = 70°, 
for photon energies E = 0.5, 3 keV. The solid curves show the results of numerical integration of the transfer equations for the Stokes 
parameters, while the dotted curves are calculated using the approximate method of § 15.21 



6 DISCUSSION 

We have presented a new method for incorporating partial conversion of photon modes due to vacuum polarization into 
fully-ionized, self-consistent atmosphere models of magnetized NSs. This method takes into account the nontrivial probability 
of photon mode conversion at the vacuum resonance. While recent works have clearly identified the important effects of the 
vacuum resonance and related mode conversion in determining the atmosphere radiation spectrum and polarization (Lai & 
Ho 2003a), so far the implementation of these effects in self-consistent atmosphere models, for technical reasons, has only 
considered two extreme limits: complete mode conversion and no mode conversion (Ho & Lai 2003,2004). With a direct, semi- 
explict Runga-Kutta integration of the radiative transfer equations for the photon modes (as opposed to the forward-backward 
substitution procedure of the Feautrier method) and the use of an accurate mode conversion formula for each photon, our 
new atmosphere code displays excellent stability with respect to grid resolution. Moreover, integration of the full transfer 
equations for the radiation Stokes parameters shows that our treatment of partial conversion is accurate. As expected, the 
partial conversion solution is intermediate between the extreme cases of complete and no conversion considered previously. An 
accurate treatment of vacuum polarization is a critical step toward interpreting the spectra and predicting the polarization 
signals of magnetic NSs. 

With our new atmosphere code, we have constructed a large number of atmosphere models for various magnetic field 
strengths, ranging from 10^"^ G to 5 x 10^'* G, for both H and He compositions. In agreement with previous, approximate 
calculations (Ho & Lai 2003,2004), we find that for B ^ 7 x 10^^ G the vacuum resonance affects the atmosphere spectra 
(i.e., the hard spectral tail due to nongrey opacities tends to be suppressed and the spectral line EW is reduced), and the 
effects become more significant as the magnetic field strength is increased. For B <, 7 x 10^'' G, the effect of the vacuum 
resonance on the spectra is smaller and becomes neglig ible for B ^ 4 X lO" G. However, even for such "low" field strengths, 
vacuum resonace has a significant effect on the observed X-ray polarizations (see Lai & Ho 2003b). Our new calculations 
presented in this paper are particularly important for the "intermediate" field regime (4 x lO'^^ ^ B ^ 2 x 10^'* G), where 
previous approximate treatments are inadequate. For the first time, we are able to accurately determine the B-dependence of 
the structure, spectra and polarization signals of ionized NS atmospheres. 

Since the most time-consuming and difficult part of atmosphere modeling involves finding the atmosphere temperature 
profile that satisfies the condition of radiative equilibrium, for the convenience of the astrophysics community, we have 
presented fitting formulae for the temperature profiles of various atmosphere models (see >I4.1I| . With these analytic expressions, 
it is relatively straightforward (using the procedure outlined in tl3.1.4|l to calculate various properties of the emergent radiation. 
These analytic profiles will also be useful for comparison with future theoretical atmosphere models. 

We note that the models presented in this paper have several limitations: (1) The models assume that the magnetic 
field lies along the surface normal. While a more general magnetic field inclination does not change the main results of our 
paper (e.g. the effect of vacuum resonance and the dependence of the atmosphere spectra on B), to confront observations (see 
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below), synthetic spectra must be constructed using realistic magnetic field and surface temperature distributions, adding up 
contributions over the entire NS surface. Such calculations are necessarily model-dependent, but they are needed for proper 
interpretation of observations.^ 

(2) At high density, the radiative transfer equation breaks down, due to the dense plasma effect. At large optical depth, 
the photon polarization develops a non- negligible longitudinal component and the index of refraction deviates significantly 
from unity. This occurs when the plasma frequency of the medium exceeds the photon frequency. To date, no detail ed studies 
of the transfer of radiation in dense plasmas have been performed, though it has been treated in an ad-hoc way bv lHo et alJ 
i2003l) . Nevertheless, this effect is important for treating thermal radation in the optical band, and for magnetars can affect 
the emission spectrum in the soft X-ray 1 keV) . (3) The assumption of fully ionized atmospheres may not be valid for cool 
NSs (such as the dim isolated NSs) or even the higher temperature AXPs and SGRs. Nevertheless, we expect features due to 
bound-bound and bound-free transitions of neutral species to be suppressed in the same manner as the ion cyclotron feature 
(see Ho et al. 2003; Potekhin et al. 2004,2005 for recent works on partially ionized magnetic atmosphere models). 



6.1 Implications for observations of isolated neutron stars 

As mentioned in §1, recent observations by Chandra and XMM-Newton have shown that the quiescent thermal spectra from 
AXPs and SGRs have no observable absorption features, such as the ion cyclotron line at Esi = 0.63{Z/A)Bi4 keV. As first 
pointed out by Ifo & Lai (2003), and confirmed by our more accurate calculations presented here, the inclusion of vacuum 
polarization effects provides a natural explanation for the non-detection: at _B = 5 x 10^* G, the vacuum-suppressed EW 
of the H or He cyclotron line is smaller than the current detector resolution. We expect that in the magnetar field regime, 
bound-bound and bound- free features will be similarly suppressed (see Ho et al. 2003; Potekhin et al. 2004). 

Prominent absorption lines (at 0.7 keV and 1.4 keV) have been detected from the source 1E1207. 4-520, a young neutron 
star (r ~ 2 X 1 0^ K) associated with a supernova remnant ijSanwal et al.ll2002l : iBignami et al.ll2003t IPe Luca et al.ll2004l : 
iMori et aliboosll Two viable (but tentative) identifications of these features are: (1) Ion cyclotron and atomic transitions 
of light-element (most likely He) atmospheres at B J; 10^* G (Pavlov & Bezchastnov 2005); (2) Atomic transitions of C or 
O atmospheres with B ^ 10^^ G (Mori et al. 2005). Based on our general result of line suppression in the magnetar field 
regime, we suggest that the first interpretation is unlikely to be correct, although a quantitative calculation of the atomic line 
strengths is needed to draw a firm conclusion.^" 

Absorption features have also been detected from three nearby, dim isolated NSs: RX J1308. 6-1-2127, RX J1605. 3-1-3249, 
and RX J0720.4-3125. While all three sources have similar effective temperatures (Tefr ~ 10^ K), their observed features 
occur at different e nergies and have varying equivalent widths: E ~ 0.2 — 0.3 keV with EW» 150 eV for RX J1308.6-f 2127 
jHaberl et alJboOSi) . E « 0.27 keV with EW» 40 eV for RX J0720.4-3125 faaberl et alJl2004al) . and E ^ 0.45 keV with 
EW« 80 eV for RX J1605. 3-1-3249 Jvan Kerkwiik et ahlEooi) . With a single line, it is difficuh to conclusively determine 
the true atmosphere composition. One possibility is that these features are associated with proton cyclotron resonance (with 
possible blending from atomic transitions) in a hydrogen atmosphere (Ho & Lai 2004). For RX J1308. 6-1-2127, the inferred 
magnetic field is 3 — 5 x 10^^ G, for which line suppression by vacuum resonance is ineffective. The broad width of the feature is 
consistent with that calculated by our hydrogen atmosphere models (see Fig.|5J. RX J1605. 3-1-3249 is also consistent with this 
picture: its feature corresponds to B ~ 7 x 10^^ G, and partial suppression of the line may account for its lower (by a factor 
of ~ 2) EW (see Fig. 1101 . The situation for RX J0720. 4-3125 is more complicated: the spectrum (including the line width) 
varies as a function of the rotation phase (Haberl et al. 2004a) and over a long timescale (a few years) (see Haberl et al. 2006). 
If its absorption feature is a proton cyclotron line, then the inferred magnetic field is too low for vacuum polarization effects 
to alter the line strength. Its small EW (40 eV) could arise if the line-emitting region (where B < 10^* G) is a small fraction 
of the NS surface; most of the surface would have B > 10^* G, requiring a highly non-dipolar surface field. Alternatively, if 
the atmosphere of RX J0720. 4-3125 is composed of helium, the required field strength is _B ~ 9 x 10^"^ G, strong enough for 
the vacuum effects to reduce the line width. 

Several recent papers have identified similar absorptio n features in other dim isolated NSs, though these observations may 
require independent confirmation and/or better statistics. iHaberl et aP l)2004bh report sp ectral featu r es for RX J0806.4-4123 
{E 0.4 - 0.46 keV; EW« 33 - 56) and RX J0420.0-5022 {E 0.3 keV; EW« 45), while IZane et alJ (|200R| ) report a spectral 
feature for RX J2143. 7-1-0654 {E « 0.75 keV; EW« 27). These features are similar to those described in more detail above, 
and suffer from the same difficulties of identification. Needless to say, since the magnetic field strengths of these NSs likely lie 



^ Synthetic spectra from NSs with realistic temperature distributions may provide a solution to the problem of the XDINS optical excess 
discussed in e.£.. |Pons et a l. 2002; Burwitz et al. 2003). 

The energy spacing of the lines at 0.7, 1.4, 2.1, and 2.8 keV is strongly suggestive of cyclotron harmonics. However, the analysis 
ofj^Mori ct al. (2005) indicates that, in all likelihood, the features at 2.1 and 2.8 keV are not real. Furthermore, a pure cyclotron line 
interpretation of these features is problematic because the strengths of the harmonics (either for electrons at 10^" G or protons at 10^'* G) 
are expected to be negligible. 
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in the range 5 x 10^ — 10^ G, for which accurate treatment of the vacuum resonance effect is crucial, the atmosphere models 
developed in this paper will be particularly useful, especially when combined with detailed modeling of (phase-dependent) 
synthetic spectra and (energy- dependent) lightcurves. 

6.2 Implications for future works 

ft is clear that further theoretical modeling of NS surface emission is needed to confront observations. Our discussion above 
(§6.1) also suggests that even accurate theoretical models and high-quality data may still be inadequate to break some of the 
degeneracy (e.g., magnetic field strength and geomery, atmosphere composition, surface temperature distribution) inherent 
in the problem. In this regard. X-ray polarimetry is highly desirable. Our calculations in ^l^show that polarization signals 
are complementary to X-ray spectra. For example, the polarization signals of magnetars and NSs with "ordinary" field 
strengths are qualitatively different. It is possible for a NS with a "boring" spectrum and lightcurve to generate an interesting 
polarization signature! 
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